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ABSTRACT 

This  study  examines  the  performance  of  an  axisymmetric  nozzle  which  was  designed 
to  produce  uniform,  parallel  flow  with  a  nominal  Mach  number  of  8.  A  free-piston-driven 
shock  tube  was  used  to  supply  the  nozzle  with  high-temperature,  high-pressure  test  gas. 
Performance  was  assessed  by  measuring  Pitot  pressures  across  the  exit  plane  of  the  nozzle 
and,  over  the  range  of  operating  conditions  examined,  the  nozzle  produced  satisfactorily  test 
flows.  However,  there  were  flow  disturbances  that  persisted  for  significant  times  after  flow 
initiation. 

The  detailed  starting  process  of  the  nozzle  was  also  investigated  by  performing  numeri¬ 
cal  simulations  at  several  nominal  test  conditions.  The  classical  description  of  the  starting 
process,  based  on  a  quasi-one-dimensional  model,  provided  a  reasonable  approximation  and 
was  used  to  demonstrate  that  the  starting  process  could  consume  a  significant  fraction  of 
the  otherwise  usable  test  gas.  This  was  especially  important  at  high  operating  enthalpies 
where  nozzle  supply  conditions  were  maintained  for  shorter  times.  Multidimensional  simu¬ 
lations  illustrated  a  mechanism  by  which  the  starting  process  in  the  actual  nozzle  could  take 
longer  than  that  predicted  by  the  quasi-one-dimensional  analysis.  However,  the  cause  of  the 
persistent  disturbances  observed  in  the  experimental  calibration  was  not  identified. 


Research  was  conducted  under  NASA  Contract  No.  NAS1-18605  while  the  author  was  in  residence 
at  the  Institute  for  Computer  Applications  in  Science  and  Engineering  (ICASE),  NASA  Langley  Research 
Center,  Hampton,  VA  23665. 
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1  INTRODUCTION 


In  recent  years  there  has  been  a  renewal  of  interest  in  pulse-type  aerodynamic  test  facilities 
as  a  way  of  providing  experimental  data  for  flight  Mach  numbers  above  8.  In  particular, 
the  free-piston-driven  reflected-shock  tunnel  has  emerged  as  a  versatile  facility  covering  the 
hypervelocity  range  10  <  Mfught  <  25  with  sufficient  test  section  density  to  be  useful  in  com¬ 
bustion  studies.  See  Stalker  [1]  for  a  review  of  pulse  facilities  and  their  use  in  experimental 
hypervelocity  aerodynamics. 

When  designing  experiments  for  these  facilities,  the  duration  of  the  quasi-steady  test 
flow  is  often  a  limiting  factor.  Typical  test  times  range  from  ^  to  3  millseconds.  From  the 
experimenter’s  point  of  view,  the  model  size  must  be  restricted  such  that  the  important  flow 
features  are  allowed  to  reach  steady  state  in  the  available  test  time  (see  e.g.  [2]).  From  the 
facility  designer’s  point  of  view,  it  is  important  to  make  the  steady  flow  duration  as  long  as 
possible.  Therefore,  the  facility  nozzle  must  not  only  provide  uniform  and  parallel  flow  to 
the  test  region,  it  should  also  start  quickly  and  reach  steady  state  while  consuming  as  little 
test  gas  as  possible.  These  constraints  become  more  severe  at  the  higher  enthalpies  where 
the  facility  supplies  a  smaller  amount  of  test  gas  for  a  shorter  time. 

This  paper  consists  of  two  major  parts.  The  early  sections  (2  -  4)  consider  the  design 
and  calibration  of  an  antisymmetric  nozzle  for  the  T4  shock  tunnel  facility  [3]  located  at  the 
University  of  Queensland,  Australia.  The  nozzle,  intended  to  supply  test  gas  to  a  model 
scramjet  combustor,  was  designed  using  the  simple  flow  decoupling  approach  described  in 
[4].  The  fabricated  nozzle  was  1.86m  long,  had  a  15.24mm  diameter  throat  and  a  388 mm 
diameter  exit.  Its  design  Mach  number  was  8.  The  nozzle  was  calibrated  by  measuring 
Pitot  pressures  near  the  exit  plane  for  a  number  of  operating  conditions.  In  the  calibration 
experiment,  it  was  found  that  the  nozzle  flow  was  taking  much  longer  to  reach  steady  state 
than  initially  expected. 


The  later  sections  of  this  paper  (5-7)  describe  attempts  to  numerically  simulate  the 
transient  starting  process  of  the  nozzle.  The  aim  of  the  simulations  was  to  identify  the  mech¬ 
anism  which  caused  the  starting  delay  observed  in  the  experiment.  A  quasi-one-dimensional 
method  was  used  to  investigate  the  transient  wave  structure  and  the  nozzle  starting  time 
over  a  number  of  tunnel  operating  conditions  while  an  axisymmetric  Navier-Stokes  code  was 
used  to  illustrate  the  multidimensional  and  viscous  aspects  of  the  starting  process. 
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2  REFLECTED-SHOCK  TUNNEL  OPERATION 


The  principal  features  of  a  free-piston  driven  shock  tunnel,  along  with  an  approximate  wave 
diagram,  are  shown  in  Fig.  1.  The  driver  tube,  which  initially  contains  low  pressure  helium 
downstream  of  the  piston,  and  the  shock  tube  which  contains  the  test  gas,  are  separated 
by  the  primary  diaphragm.  This  diaphragm  was  typically  composed  of  two  sheets  of  mild 
steel.  Attached  to  the  downstream  end  of  the  shock  tube  is  the  facility  nozzle  whose  throat 
is  significantly  smaller  than  the  inside  diameter  of  the  shock  tube.  The  subsonic  portion 
of  the  nozzle  effectively  closes  the  downstream  end  of  the  shock  tube  and  forms  the  shock 
reflection  region.  The  supersonic  portion  of  the  nozzle  empties  directly  into  a  test  section 
and  dump  tank  which  is  evacuated  to  an  initial  pressure  of  approximately  30  Pa.  The  test 
gas  is  retained  in  the  shock  tube  by  a  thin  mylar  diaphragm. 

The  first  stage  of  operation  is  the  launch  of  the  piston  and  its  acceleration  along  the 
compression  tube.  The  driving  force  is  supplied  by  compressed  air  from  a  reservoir.  The 
helium  in  front  of  the  piston  is  compressed  and  eventually  reaches  the  same  pressure  as  the 
compressed  air  driving  the  piston.  Because  of  its  large  kinetic  energy,  the  piston  overshoots 
and  further  compresses  the  helium.  The  primary  diaphragm  subsequently  bursts  (at  a  pres¬ 
sure  5 6.6MPa  for  a  4mm  thick  diaphragm).  At  this  point,  the  helium  has  been  compressed 
to  1  / 60th  of  its  initial  volume  and  is  contained  in  approximately  0.5m  of  the  compression 
tube.  This  part  of  the  compression  tube  then  becomes  a  constant  pressure  driver  section  for 
the  shock  tube  [5]. 

After  rupture  of  the  primary  diaphragm,  the  helium  expands  into  the  shock  tube  and 
compresses  the  test  gas  before  it.  Unlike  the  approximately  isentropic  compression  of  the 
helium,  the  processing  of  the  test  gas  is  performed  by  a  strong  normal  shock.  This  primary 
shock  wave  travels  the  length  of  the  shock  tube,  reflects  from  the  closed  end,  and  brings  the 
test  gas  to  rest  in  the  nozzle  supply  region.  Operation  in  this  manner  is  called  tailored  [6] 
and  is  shown  in  the  wave  diagram  (Fig.  1(b))  by  the  contact  surface  coming  to  rest  when 
intercepted  by  the  reflected  shock.  The  shock-compressed  test  gas  is  temporarily  contained 
in  the  nozzle  supply  region  which  has  a  length  of  0.25  —  0.5m.  Ideally  the  nozzle  supply 
conditions,  characterized  by  the  total  enthalpy  H,  and  pressure  P,,  are  maintained  as  the 
reflected  shock  continues  upstream  though  the  driver  gas. 

Aside  from  supply  pressure  fluctuations  after  shock  reflection,  the  useful  test  time  is 
terminated  by  contamination  of  the  test  gas  by  driver  gas.  The  mechanism  for  this  con¬ 
tamination  is  the  bifurcation  of  the  reflected  shock  into  two  oblique  shocks  near  the  tube 
wall  [7]  [8].  While  the  strong  normal  shock  brings  the  driver  gas  to  rest  m  the  centre  of 
the  tube,  the  weaker  oblique  shocks  allow  gas  to  jet  along  the  tube  walls  and  prematurely 
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reach  the  nozzle  supply  region.  In  an  effort  to  delay  the  arrival  of  this  jet,  the  shock  tube 
is  usually  operated  in  an  undertailored  mode  where  the  reflected  shock  accelerates  into  the 
driver  gas  and  an  expansion  propagates  into  the  nozzle  supply  region.  This  increases  the 
distance  between  the  driver-gas/test-gas  interface  and  the  end  wall  of  the  shock  tube  but 
also  results  in  an  unavoidable  drop  in  P,.  The  approximate  time  at  which  contamination  of 
the  test  gas  is  expected  to  occur  at  the  nozzle  throat  decreases  with  increasing  stagnation 
enthalpy  (see  [1]). 

Although  the  ideal  operation  of  the  facility  provides  constant  pressure  driver  conditions 
to  the  shock  tube,  the  effect  of  finite  driver  size  and  non-optimal  piston  dynamics  on  the 
relaxation  of  P,  is  noticeable.  These  effects  result  in  the  long  term  relaxation  of  P,  and 
may  also  contribute  to  the  early  decay  of  Fs.  The  net  result  is  a  typical  nozzle  supply 
pressure  history  as  shown  in  Fig.  2.  The  transducer  used  to  obtain  this  trace  was  located 
approximately  8  centimetres  upstream  of  the  closed  end  of  the  shock  tube.  Hence,  the 
passage  of  the  primary  and  reflected  shocks  are  shown  as  distinct  events.  Because  of  the 
location  and  limited  response  time  for  the  transducer,  the  peak  reflection  pressure  was  not 
recorded.  Past  the  maximum  value,  P,  continued  to  decay  due  to  the  combined  effects  of 
undertailored  operation  and  driver  dynamics.  For  pure  helium  driver  gas  and  operation 
considered  here,  this  decay  was  typically  25  —  30%  during  a  nominal  0.5ms  test  time. 

Upon  shock  reflection,  the  light  secondary  diaphragm  bursts  and  some  of  the  test  gas 
following  the  primary  shock  expands  through  the  nozzle  throat  into  the  divergent  part  of  the 
nozzle.  From  the  point  of  view  of  the  nozzle,  the  shock  tube  is  now  a  reservoir  of  stagnant, 
high-temperature,  high-pressure  test  gas. 

The  subsequent  starting  of  the  nozzle  has  been  examined  experimentally  and  analytically 
by  Amann  [9]  and  Smith  [10]  [11].  As  the  primary  shock  travels  down  the  nozzle  it  acceler¬ 
ates  the  very  low  pressure  gas  already  in  the  nozzle.  However,  the  primary  shock  decelerates 
because  of  the  diverging  nozzle  walls  [12].  Figure  3,  based  on  the  quasi-one-dimensional 
model  in  [11]  illustrates  this  situation.  Test  gas  which  accelerates  through  the  nozzle  throat 
following  the  primary  shock,  expands  to  a  very  high  Mach  number  and  is  suddenly  decel¬ 
erated  when  it  encounters  the  slower  primary  shock  structure.  An  upstream-facing  shock 
is  thus  formed  and  is  swept  downstream  through  the  nozzle.  Between  the  upstream-facing 
shock  and  the  steady  expansion  being  established  near  the  throat  there  is  an  unsteady  ex¬ 
pansion,  the  upstream  head  of  which  is  being  swept  downstream  with  velocity  u  —  a.  When 
the  unsteady  expansion  is  the  last  wave  to  be  expelled,  the  starting  time  for  the  nozzle 
can  be  approximated  from  the  steady  state  conditions.  For  higher  initial  pressures  in  the 
nozzle/test-section,  the  transit  time  for  the  upstream-facing  shock  will  determine  the  noz- 
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zle  starting  time  [10]  [11].  In  the  conical  and  wedge-shaped  nozzles  studied  by  Smith,  the 
initial  pressure  could  be  significantly  higher  than  the  steady-state  static  pressure  without 
prolonging  the  starting  processes.  However,  for  the  Mach  8  contoured  nozzle  studied  here, 
an  initial  pressure  of  a  tenth  of  the  steady-state  static  pressure  is  sufficient  to  delay  the 
quasi-one- dimensional  starting  process. 

The  quasi-one-dimensional  model  was  considered  to  be  adequate  for  most  practical  situ¬ 
ations  where  the  divergence  angle  of  the  nozzle  is  small  and  the  design  Mach  number  of  the 
nozzle  is  not  too  large.  Britan  and  Vasil’ev  [13],  [14]  have  studied  two-dimensional  situations 
in  which  the  nozzle  divergence  angle  is  large  and  where  multidimensional  effects  upstream 
of  the  throat  result  in  high  peak  properties  at  the  throat  (on  a  time  scale  much  smaller  that 
the  nominal  test  time).  This  study  will  not  consider  such  effects  in  the  subsonic  part  of  the 
nozzle  but  will  consider  the  situation  where  the  nozzle  is  relatively  long  (when  compared 
to  the  equivalent  conical  nozzle  with  the  same  initial  divergence),  has  a  significantly  ihick 
boundary  layer,  and  is  operated  near  conditions  where  the  upstream-facing  shock  limits  the 
starting  time. 


3  NOZZLE  DESIGN 

Once  the  test  gas  is  brought  to  rest  at  the  end  of  the  shock  tube,  it  is  expanded  to  the  test 
section  conditions  in  two  conceptually  separate  stages.  First,  the  gas  is  expanded  through  a 
conical  section  to  produce  a  uniform  source  flow  and  then  it  is  straightened  by  the  contoured 
section  of  the  nozzle  wall  such  that  the  flow  into  the  test  section  is  approximately  uniform 
and  parallel. 

The  design  point  chosen  for  the  Mach  8  nozzle  had  a  nozzle  supply  pressure  P,  = 
30.4MPa,  a  supply  enthalpy  H,  =  16.2 MJ/kg  and  temperature  T,  =  8000/f.  The  throat 
diameter  d,  —  15.24mm  was  chosen  after  some  iteration.  A  conceptual  view  of  the  entire 
nozzle  is  shown  in  Fig.  4.  For  computational  convenience,  we  set  the  origin  of  the  axial 
coordinate  to  the  start  of  the  conical  expansion  (immediately  after  the  throat). 

For  the  subsonic  portion  of  the  nozzle,  most  of  the  end  of  the  shock-tube  had  been  kept 
flat  and  normal  to  the  tube  axis  so  that  the  shock  reflection  was  well  behaved.  This  was  not 
expected  to  cause  any  problem  as  the  area  ratio  of  the  nozzle  was  large  ( Aextt/A *  ~  667) 
[15].  The  subsonic  lead-in  to  the  throat  needed  to  be  smooth  in  order  to  avoid  separation 
and  associated  oblique  shocks  at  the  beginning  of  the  conical  expansion  [9].  The  throat  was 
a  constant  diameter  section  of  length  d»  and  was  manufactured  from  a  Beryllium-Copper 
alloy  (BIO  or  B25).  In  reference  [16]  it  was  indicated  that  such  throats  produce  good  quality 


4 


source  flow  when  used,  in  conical  nozzles.  Also,  it  was  found  that  these  throats  resisted 
erosion  better  than  those  with  a  shorter  constant-diameter  section  or  those  manufactured 
from  grade  316  stainless  steel.  In  operation,  the  Be-Cu  throats  increased  in  diameter  by 
approximately  0.5mm  over  six  shots  when  the  tunnel  was  operated  at  H,  =  lOMJ/kg  and 
P.  =  50  -  70 MPa. 

The  design  of  the  initial  conical  expansion  (region  2,  Fig.  4)  involved  the  analysis  of 
the  flow  as  a  steady,  quasi-one-dimensional  chemically  reacting  flow.  Flow  properties  were 
determined  with  the  program  NENZF  [17]  in  which  gas  in  the  nozzle  supply  region  (region 
1)  upstream  of  the  throat  is  assumed  to  be  stagnant  and  in  chemical  equilibrium.  A  chemical 
equilibrium  model  was  also  used  in  the  throat  but  a  finite-rate  chemistry  model  was  used 
downstream  of  the  throat.  The  conical  expansion  had  a  cone  half-angle  of  14°  and  extended 
to  x  =  0.168m.  At  this  point,  the  reactions  were  assumed  frozen  and  the  gas  was  modelled 
as  a  perfect  gas  with  a  suitably  chosen  ratio  of  specific  heats,  7.  The  flow  conditions  at  the 
end  of  the  conical  expansion  were  M \  =  it/ a  —  4.245  where  a  is  an  estimate  of  the  speed 
of  sound.  The  value  7  =  1.38  was  chosen  by  considering  equivalent  chemically  reacting  and 
perfect  gas  flows  through  the  same  quasi-one-dimensional  expansion. 

The  Method-of- Characteristics  (MOC)  was  then  used  to  compute  the  source  flow  (region 
3),  the  flow  correcting  section  (region  4)  and  the  uniform  flow  region  (region  5)  in  a  semi- 
automated  way  with  the  program  described  in  [18].  The  inlet  boundary  was  specified  as 
uniform  source  flow  at  12  points  along  the  arc  AB.  The  Mach  number  on  the  axis  at  point 
C  was  computed  to  be  Me  =  8.02  and,  beyond  this  point,  was  held  constant.  The  solution 
then  stepped  downstream  along  the  axis  and  computed  the  required  flow  in  region  4  by 
proceeding  upstream  along  characteristics  of  the  same  family  as  CA.  Once  the  characteristic 
mesh  (Fig.  5)  was  generated,  a  streamline  was  interpolated  through  the  mesh,  starting  at 
point  A  and  finishing  where  it  intersected  the  characteristic  CD  (which  forms  the  upstream 
boundary  of  the  uniform  test  flow  region).  The  data  points  on  the  interpolated  streamline 
were  then  fitted  with  a  cubic-spline  function  [19]  having  eight  knots  and  a  starting  slope 
specified  to  match  the  initial  conical  expansion.  Table  1  displays  the  coordinates  of  the 
knots  for  the  spline  function. 

The  design  process  described  above  produces  a  nozzle  with  the  smallest  inviscid  design 
length  Ld  =  2.187m.  However,  a  smoother  transition  may  be  obtained  by  adding  a  section 
with  a  smooth  variation  of  axial  Mach  number  between  the  end  of  the  source  flow  region 
(C)  and  the  start  of  the  uniform  test  flow  region.  The  length  of  the  fabricated  nozzle  was 
later  reduced  by  truncating  the  design  contour  at  Lt  =  1.8m  where  the  characteristic  CD 
intersected  the  estimated  boundary  layer  edge.  Because  of  the  high  exit  Mach  number,  Lt 
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j 

Xj 

(m) 

rj 

(m) 

comment 

1 

0.16800 

0.04951 

dr/dx  =  0.2493 

2 

0.45648 

0.10584 

3 

0.74488 

0.14208 

4 

1.03336 

0.16552 

5 

1.32176 

0.18104 

6 

1.61024 

0.19040 

7 

1.89872 

0.19520 

8 

2.18712 

0.19688 

natural  end  condition 

Table  1:  Knots  for  the  cubic  interpolating  spline, 
is  significantly  smaller  than  Ld- 

The  final  design  is  shown  in  cross-section  in  Fig.  6.  Note  that  the  converging  (subsonic) 
region  closes  off  the  shock  tube  to  form  the  shock  reflection  region.  This  component  of  the 
nozzle  is  manufactured  from  high-tensile  steel  and  extends  into  the  diverging  section  of  the 
nozzle  to  the  axial  location  x  =  0.18m.  Most  of  the  contoured  wall  is  manufactured  from 
glass-reinforced  plastic  and  is  attached  to  the  steel  section  with  an  aluminium  flange. 


4  EXPERIMENTAL  CALIBRATION 

The  performance  of  the  nozzle  was  evaluated  by  measuring  the  Pitot  pressure  Pjntot  at  planes 
normal  to  the  nozzle  axis  and  located  at  x  ~  1.8m  and  x  ~  2.0m.  Each  Pitot  probe  was  fitted 
with  a  PCB-112  piezo-electric  pressure  transducer  which  measured  the  stagnation  pressure 
behind  a  detached  shock  that  formed  over  the  upstream  face  of  the  probe.  Several  probes 
were  mounted  in  a  rake  and  a  number  of  shots  of  the  shock  tunnel  were  made  at  each  nominal 
operating  condition.  Refer  to  [20]  for  data  recording  and  analysis  details  and  a  complete  set 
of  data. 

Figure  7  shows  the  unfiltered  histories  of  both  nozzle  supply  pressure  and  Pitot  pressure 
on  the  nozzle  axis  at  x  ~  1.8m.  These  shots  were  typical  of  the  ~  30  shots  that  were  made 
for  this  survey.  Note  that  the  time  origin  is  arbitrary  as  the  recording  device  was  triggered 
from  the  supply  pressure  signal  and  had  an  arbitrary  pretrigger  delay.  All  of  the  traces  in 
Fig.  7  show  the  impulsive  start  and  subsequent  decay  as  discussed  in  Section  2.  The  time 
delay  from  shock  reflection  to  the  initial  rise  in  the  Pitot  trace  is  the  time  for  the  primary 
shock  to  propagate  through  the  nozzle.  The  relatively  slow  initial  rise  of  the  Pitot  traces  is 
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due  to  the  filling  of  the  probe  cavity. 

Filtered  Pitot  pressure  traces  at  a  number  of  radial  locations  are  shown  in  Fig.  8.  Here, 
the  time  base  has  been  shifted  such  that  the  shock  reflection  occurred  at  t  =  0.  The  filtering 
was  done  taking  a  moving  average  of  the  trace  over  a  0.05ms  window.  Note  that,  after 
reaching  peak  values,  the  traces  near  the  nozzle  axis  show  a  rapid  decay  down  to  the  levels 
of  the  outer  traces.  Note  also  the  fall-off  in  average  pressure  and  the  decline  in  quality  for 
the  traces  at  r  >  140mm.  This  location  appears  to  be  within  the  boundary  layer  formed 
along  the  nozzle  wall. 

Although  the  absolute  pressures  are  changing  during  the  test  flow  period,  a  quasi-steady 
estimate  can  be  obtained  by  considering  Pitot  traces  normalized  by  a  time-shifted  supply 
pressure  trace 

P-norm  =  Ppitot{t) / P.{t  -  St)  ,  (1) 

where  St  is  the  nominal  transit  time  of  a  fluid  particle  from  the  shock  reflection  region  to 
the  probe  position.  Normalized  Pitot  pressure  histories  for  shots  1405  and  1406  are  shown 
in  Fig.  9  where  St  —  0.50ms.  The  extreme  values  in  the  early  part  of  the  traces  ( t  <  0.5ms) 
are  not  significant  as  the  nominal  time  used  in  equation  1  is  less  than  the  time  taken  for  the 
starting  shock  structure  to  traverse  the  nozzle. 

The  traces  at  r  =  112, 120mm  illustrate  the  classical  starting  process  in  which  the  shock 
structure  appears  first  followed  by  the  unsteady  expansion  fan  (the  dip  at  t  —  0.6ms  and 
subsequent  rise  to  t  —  1.0ms)  and  then  steady  conditions  [t  >  1.2ms).  Toward  the  nozzle  axis 
the  flow  establishment  is  qualitatively  different.  Starting  with  the  passage  of  the  unsteady 
expansion  ( t  ~  0.75ms),  there  is  a  large  disturbance  in  the  flow  which  is  evident  out  to  a 
radius  r  =  90mm.  Such  behaviour  was  obvious  in  most  (but  not  all  shots)  and  tended  to 
be  less  significant  at  the  highest  enthalpy  condition  ( H ,  ~  1 6MJ/kg).  For  H,  >  10 MJ/kg, 
the  disturbed  flow  consumed  a  significant  fraction  of  the  available  test  time. 

Figure  10  shows  the  tuning  of  the  events  at  the  exit  plane  of  the  nozzle  for  shots  with 
a  nominal  supply  pressure  P,  ~  50 MPa.  Predictions  for  the  arrival  time  of  the  primary 
shock  (line  A)  and  the  time  to  settle  within  10%  of  steady  state  (line  B)  were  computed 
using  the  quasi-one-dimensional  model  and  the  code  described  in  Section  6.  The  initial 
pressure  within  the  nozzle  was  set  to  P,mt  =  33Pa  which  is  typical  of  conditions  used  in  the 
T4  facility.  Unfortunately  precise  values  of  Pintl  were  not  recorded.  Approximate  times  for 
arrival  of  contaminated  gas  were  taken  from  [1]  and  adjusted  for  the  length  of  this  nozzle. 

The  arrival  time  of  the  incident  shock  was  well  approximated  by  the  quasi-one-dimensional 
model  except  at  the  H,  ~  6.6 MJ/kg  condition.  Here,  either  transient  supply  conditions  (see 
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Section  5.2)  or  higher  initial  pressures  in  the  test  section  may  account  for  the  differences. 

The  settling  time  estimates  (line  B,  Fig.  10)  were  relatively  insensitive  to  either  of  these 
influences.  In  the  experimental  traces,  the  passage  of  the  upstream  head  of  the  expansion 
fan  was  identified  as  the  point  (after  the  shock  structure  passes)  where  the  normalized  Pitot 
pressure  rose  to  the  steady  level.  (Note  that  the  trace  may  have  exceeded  the  steady-state 
value  after  this  time.)  These  times  agree  qualitatively  with  the  predictions  of  the  settling 
times  (line  B)  but  were  all  consistently  delayed.  As  will  be  seen  in  Section  7.1,  this  may  be  a 
multidimensional  effect  associated  with  the  distortion  of  the  upstream-facing  shock  and  the 
unsteady  expansion. 

The  experimental  settling  time  seemed  to  occur  approximately  0.3ms  after  the  passage 
of  the  unsteady  expansion.  At  the  highest  enthalpy  ( H a  ~  16 A/ J/ kg)  this  settling  time 
coincided  with  the  estimated  arrival  time  of  the  driver  gas  contamination  however,  the  con¬ 
tamination  time  estimates  were  considered  conservative  and  the  flow  disturbances  were  rel¬ 
atively  small  at  this  condition.  Hence,  some  useful  (albeit  short)  test  time  was  expected  to 
be  available. 

Once  a  test  time  was  identified  (by  starting  after  the  initial  flow  perturbations  and  then 
terminating  before  the  expected  time  of  contamination),  a  quasi-steady  value  was  computed 
as  the  mean  value  of  the  filtered  and  normalized  trace  over  the  test  period.  These  mean 
values,  together  with  standard  deviations  computed  over  the  same  period,  are  shown  in  Fig. 

11.  The  profiles  are  grouped  with  increasing  supply  enthalpy  to  the  right  and  increasing 
supply  pressure  to  the  top.  Other  than  some  spurious  data  points  near  the  nozzle  axis,  the 
profiles  appear  to  be  satisfactory  (for  shock  tunnel  work).  If  the  nozzle  contour  was  not 
performing  well  then,  disturbances  would  be  expected  to  focus  near  ihe  nozzle  axis.  The 
stray  data  points  may  indicate  the  existence  of  such  problems. 

5  FLOW  MODELLING 

The  goal  of  the  numerical  simulations  was  to  provide  some  insight  to  the  pertormance  of  the 
nozzle  in  terms  of  its  time  to  start  and  its  steady-state  flow  uniformity.  The  causes  of  the 
late-time  disturbances  seen  in  the  Pitot  pressure  traces  were  of  particular  interest. 

To  achieve  this  goal,  time-dependent  simulations  were  performed  with  a  quasi-one-dimensional 
(Q1D)  code  and  with  an  axisymmetric  full-Navier-Stokes  (FNS)  code  [21],  The  former  pro¬ 
vided  a  computationally  inexpensive  way  to  explore  part  of  the  parameter  space  for  the 
nozzle  operating  conditions  while  the  FNS  code  was  used  to  obtain  detailed  information  on 
the  viscous  and  multidimensional  effects  at  a  single  operating  condition. 
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In  both  approaches  the  nozzle  flow  was  simplified  as  follows.  The  gas  in  the  shock- 
reflection  /  nozzle- .supply  region  of  the  shock  tube  was  assumed  to  be  stagnant  and  in 
chemical  eq'j’rjrium.  However,  supply  conditions  were  allowed  to  change  with  time  as 
governed  by  the  supply  pressure  trace  (Fig.  2).  Flow  into  the  nozzle  was  assumed  inviscid, 
quasi-steady  and  also  in  chemical  equilibrium.  There  were  no  transient  or  shock-focusing 
effects  [13]  considered  here.  Starting  at  the  throat,  the  gas  was  assumed  to  be  a  non-reacting 
mixture  of  nitrogen  and  oxygen  with  mass  fractions  0.7686  and  0.2314,  respectively.  The 
thermodynamic  model  [22]  was  based  on  weighting  the  polynomial  fits  for  the  individual 
species.  The  particular  constants  and  expressions  used  for  the  Q1D  simulations  are  given  in 
Appendix  A.  A  further  simplification  was  that  only  laminar  flow’s  were  considered. 

5.1  Geometry 

Figure  12  shows  flow  geometry  considered  in  the  simulations.  The  axis  of  symmetry  is 
aligned  with  the  x-axis  and  flow  is  from  left  to  right,  entering  through  the  nozzle  throat 
at  x  —  0  and  exiting  at  the  x  =  2.187 m  boundary.  The  inflow  boundary  was  considered 
to  have  uniform,  supersonic  flow  conditions  as  specified  in  Section  5.2,  while  the  outflow 
boundary  condition  was  implemented  as  a  zero-order  extrapolation.  The  x  =  1.8m  plane 
was  considered  far  enough  upstream  of  the  outflow  boundary  to  be  unaffected  by  the  form 
of  the  outflow  boundary  condition,  even  in  the  presence  of  a  boundary  layer. 

The  radial  position  of  the  wall  starts  with  a  value  rwau  =  0.00762m  at  x  =  0,  varied 
linearly  with  slope  dr/dx  —  0.2493  for  0  <  x  <  0.168m  and  was  defined  as  a  cubic-spline 
function  (see  Table  1)  for  0.168m  <  x  <  2.187m.  Note  that  the  contoured  wall  of  the 
fabricated  nozzle  was  truncated  at  x  =  1.8m.  For  the  Q1D  calculations,  the  nozzle  wall  was 
considered  to  be  an  inviscid  boundary  while,  for  the  FNS  calculations,  it  was  considered  to 
be  either  the  same  or  a  no-slip  constant-temperature  ( Twau  —  300/f )  boundary.  A  constant- 
temperature  boundary  was  considered  appropriate  because  the  surface  temperature  of  the 
wall  was  expected  to  rise  no  more  than  a  few  degrees  (say,  <  10 A')  above  room  temperature 
for  the  duration  of  the  test  flow  [23]. 

5.2  Inflow  and  Initial  Conditions 

The  inflow  condition  at  x  =  0  was  set  to  approximate  the  throat  conditions  in  the  T4  facility 
at  the  nominal  operating  condition  with  a  4mm  mild  steel  primary  diaphragm,  a  driver  tube 
compression  ratio  X  —  60,  and  the  shock  tube  initially  filled  with  air  at  P  =  120 kPa  and 
T  =  296A.  The  flow  conditions  at  the  beginning  of  the  conical  expansion  were  obtained 
with  the  following  procedure. 
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Step  1:  The  supply  conditions  (behind  the  reflected  shock)  were  estimated  with  the 
program  ESTC  [24]  which  incorporated  an  equilibrium  chemistry  model  for  air  with  the 
species  N2,  N,  02,  O,  NO,  NO+  and  e~ .  From  the  shock  reflection  conditions,  the  test 
gas  was  allowed  to  expand  adiabatically  in  chemical  equilibrium  to  the  measured  equilibrium 
supply  pressure  [15].  From  Fig.  13,  the  equilibrium  pressure  was  estimated  to  be  P,_rey  ~ 
50 MPa  at  t  =  3.4 ms  thus  giving  an  equilibrium  stagnation  enthalpy  //,_rey  ~  10 MJ/kg 
and  temperature  T,_re/  ~  6010/ST. 

Step  2:  Using  these  reference  conditions,  flow  at  the  throat  was  estimated  using  the 
quasi-one-dimensional  code  NENZF  [17]  in  which  the  test  gas  was  again  assumed  to  be  in 
equilibrium.  Conditions  at  the  throat  consistent  with  the  N2  —  02  mixture  were:  temperature 
T —  5448 K,  pressure  P,  —  27 ,9AM Pa,  density  p .  =  17.79 kg/m3,  velocity  ut  —  1450 m/s, 
and  specific  internal  energy  e*  =  5.24 MJ/kg. 

Step  3:  A  transient  inflow  was  constructed  from  the  supply  pressure  history  by  using 
ESTC  and  NENZF  to  compute  the  quasi-steady  throat  conditions  for  a  number  of  P,  values. 
The  control  points,  labelled  a  through  d  in  Table  2,  divided  the  supply  pressure  trace  into 
three  convenient  segments,  the  first  two  of  which  were  specified  as  the  cubic 

P,  =  a3r3  -f  a2r2  +  aqr  +  a0,  (2) 

where  the  a,  and  r  are  defined  in  Table  3  and  t  <  tc.  When  t  >  tc,  the  supply  pressure 
followed  the  exponential  decay 

P,  =  Pcexp(-a4(t  -  tc)),  (3) 


where 


at=-,o9(Pd/Pc)  {4) 

\td  —  tcj 

Figure  13  shows  a  comparison  between  these  approximate  expressions  for  supply  pressure 
and  the  measured  (filtered)  supply  pressure  trace  for  shot  1405.  The  quasi-steady  supply 
temperature  was  found  to  be  well  approximated  by  the  perfect  gas  relation 

(7-l)/7 

I  (5) 


T, 


P, 


X- 


r'f, 


,  Pi-ref , 


where  the  effective  ratio  of  specific  heats  is  7  =  1.197  and  the  reference  conditions  are  the 
nominal  supply  conditions  given  above. 

Step  4:  The  conditions  at  the  nozzle  throat  were  then  related  to  the  transient  value  of 
P,  with  the  following  correlations 


P./P,  =  0.559, 
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Table  2:  Control  points  for  the  transient  Pa  expressions. 


Parameter 

ta  <t  <tb 

ib  <t  <tc 

r 

{t  ~  ta)/(tb  ~  ta ) 

( t  -  tb)/{tc  -  h) 

<*o 

Pa 

Pb 

a  i 

f (ft  -  ft) 

0 

0 

1  (ft  -  ft) 

«3 

!(ft  -  ft) 

\{Pb  -  Pc) 

Table  3:  Parameters  for  the  transient  supply  pressure  equation. 
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T./T, 

P*/ps 

u./aa 


(6) 


=  0.907, 

=  0.624, 

=  1.050, 

where  a ,  —  (317.0  *  T,)1^2 .  Other  properties  are  obtained  from  the  thermodynamic  model 
for  the  N2  —  02  mixture.  Figure  14  shr  ws  the  variation  of  these  properties  with  time.  Note 
that,  while  the  throat  pressure  and  density  drop  approximately  30%  from  t  =  0.5ms  to 
t  —  1.0ms,  the  temperature  and  velocity  drop  by  only  7.5%  and  5%,  respectively. 

Initial  conditions  in  the  nozzle  (and  facility  test  section)  were  set  to  Tinit  —  300  if ,  utm(  —  0 
and  33Pa  <  P{nit  <  150Pa.  The  lower  values  of  Pxnu  =  33Pa  were  typical  of  the  facility 
operating  conditions  but  the  higher  pressures  may  occur  in  the  presence  of  test  section  leaks 
or  when  hydrogen-fueled  scramjet  models  are  being  tested.  For  such  tests,  hydrogen  fuel  is 
injected  into  the  model  (and  test  section)  several  milliseconds  before  the  flow  through  the 
nozzle  is  started. 


6  QUASI-ONE-DIMENSIONAL  SIMULATIONS 

The  quasi-one- dimensional  simulations  were  performed  for  a  single  nozzle  supply  condition 
Ha-Tef  =  10 MJ/kg  with  P,_re/  =  50 MPa.  Initial  pressures  ranged  from  Ptn;t  =  33Pa 
to  Pinit  =  300Pa.  The  lower  values  of  P*nit  are  representative  of  actual  facility  operating 
conditions  while  the  higher  values  were  performed  for  comparison  with  the  axisymmetric 
simulations. 

The  code  used  a  simple  finite-volume  formulation  which  was  second  order  accurate  in 
both  space  and  time.  An  Osher-type  Riemann  solver  [25]  was  used  to  obtain  interface  fluxes 
and  a  SUPERBEE  limiter  [26]  was  used  to  interpolate  flow  properties.  The  implementation 
is  more  fully  described  in  Appendix  B. 

In  all  of  the  simulations  reported  here,  the  computational  cells  were  arranged  with  equal 
increments  in  x ,  and  time  stepping  was  done  with  a  CFL  number  of  approximately  0.25. 
Simulations  with  250  cells  were  used  to  obtain  the  timing  data  shown  in  Fig.  10  while 
simulations  with  500  cells  were  used  to  obtain  the  results  discussed  below.  There  was  little 
difference  (i.e.  <  1%)  between  both  sets  of  solutions  in  terms  of  both  wave  speeds  and 
steady-state  conditions. 
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6.1  Results 


Figure  15  shows  time  sequences  of  (logarithmic)  density  profiles  for  PlTllt  =  33 Pa  and  150Pa. 
Steady  inflow  conditions  were  maintained  at  the  throat.  All  of  the  major  features  of  Smith’s 
[11]  quasi-one-dimensional  model  appeared.  For  Pinit  =  33 Pa,  the  unsteady  expansion  was 
the  last  wave  to  sweep  through  the  nozzle  but,  for  Pinit  =  150 Pa,  the  unsteady  expansion 
has  merged  into  the  upstream- facing  shock  before  this  shock  reaches  x  =  1.8m.  Hence,  the 
33Pa  simulation  provided  an  example  of  the  minimum  starting  time  while  the  slower  waves 
in  the  150Po  simulation  prolonged  the  starting  process.  Note  that  PlTilt  =  150Pa  was  well 
below  the  steady-state  static  pressure  of  1667Pa  that  eventually  developed  at  the  nozzle 
exit. 

This  was  a  quite  different  behaviour  to  that  of  the  conical  and  wedge  nozzles  studied  in 
[11]  where  Pinn  could  be  significantly  larger  than  the  steady-state  static  pressure  without 
prolonging  the  starting  process.  The  key  difference  in  the  flow  geometry  is  the  (necessarily) 
long  contoured  section  of  the  nozzle  which  allows  the  unsteady  expansion  to  catch  (and 
merge  into)  the  slower  moving  shock  structure.  This  is  shown  clearly  in  Fig.  16  where  the 
wave  trajectories  (obtained  from  Fig.  15)  are  plotted  for  a  number  of  initial  pressures.  The 
trajectory  of  the  upstream  head  of  the  unsteady  expansion  is  essentially  the  same  in  all  cases 
but  the  starting  shock  structure  slows  significantly  for  higher  values  of  Pinit- 

Simulations  were  also  performed  with  transient  throat  conditions  as  described  in  Section 
5.2.  A  summary  of  times,  for  arrival  of  the  primary  shock  at  x  =  1.8m  and  the  passage 
of  the  upstream  head  of  the  unsteady  expansion  at  the  same  point  is  shown  in  Fig.  17. 
Transition  from  a  minimum  starting-time  condition  (limited  by  the  unsteady  expansion)  to 
a  condition  limited  by  the  upstream-facing  shock  occurs  in  the  range  60Pa  <  Pinit  <  120Pa. 
The  ranges  of  event  times  for  the  experimental  data  are  shown  as  the  hatched  regions. 
The  primary  shock  arrival  time  is  well  approximated  by  the  steady  inflow  simulations  (in 
preference  to  the  transient  inflow  conditions).  This  may  be  caused  by  the  omission  of  the 
very  high  peak  pressures  from  the  transient  inflow  specification  and  the  shock  focusing  effects 
examined  by  Britan  and  Vasil’ev  [13].  Both  of  these  effects  would  result  in  faster  primary 
shocks.  The  arrival  of  the  upsteam  head  of  the  unsteady  expansion  occurred  at  a  later  time 
than  that  predicted  by  either  of  the  inflow  conditions.  This  delay  may  be  attributed  to 
multidimensional  effects  as  discussed  in  Section  7.1. 

Histories  of  the  simulated  Pitot  pressure  for  both  steady  and  unsteady  throat  conditions 
are  compared  with  a  typical  experimental  measurement  in  Fig.  18.  The  time  base  for  the 
experimental  data  has  been  shifted  to  put  the  shock  reflection  (in  the  nozzle  supply  region) 
at  approximately  t  =  0.  The  Pitot  pressure  displayed  for  the  simulations  was  computed  with 
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the  approximate  relation 


Ppitot  =  0.92 pu2  ,  (7) 

which  is  suitable  for  high  Mach  number  flow.  For  the  moment,  ignore  the  differences  in 
steady  state  levels  as  these  may  be  attributed  to  boundary  layer  effects  not  included  in 
the  quasi-one-dimensional  simulations.  The  arival  times  of  primary  shock,  upstream-facing 
shock  and  expansion  fan  are  reasonably  close.  However,  the  disturbance  in  the  experimental 
data  from  t  =  0.7 ms  to  1.0ms  has  no  counterpart  in  the  quasi-one-dimensional  numerical 
simulations. 


7  AXISYMMETRIC  SIMULATIONS 

The  axisymmetric  simulations  were  performed  in  the  hope  of  identifying  multidimensional  or 
viscous  mechanisms  causing  the  persistent  disturbances  seen  in  the  Pitot  pressure  histories 
of  the  experimental  calibration  (cf.  Fig.  18b).  Two  simulations  were  performed  with  Pinit  = 
150Pa  and  the  steady  throat  conditions  described  in  Section  5.2.  One  used  an  inviscid 
wall  boundary  condition  while  the  second  used  a  no-slip,  constant  temperature  boundary 
condition. 

The  simulations  were  performed  with  an  axisymmetric  version  of  the  SPARK  code  [21] 
which  is  a  finite-difference  program  that  integrates  the  conservation  equations  for  mass, 
momentum,  energy  and  species  in  a  time-accurate  manner.  Along  the  x-axis,  the  singular 
governing  equations  were  collapsed  using  L’Hopital’s  rule.  The  option  selected  for  the  time 
integration  was  the  second-order  MacCormack  scheme  with  a  CFL  number  of  0.6.  Numerical 
damping  was  implemented  as  a  mixture  of  a  standard  4th  order  artificial  dissipation  (25%) 
and  a  total-variation- diminishing  dissipation  term  (75%)  [27]  based  on  static  pressure  only. 
The  dissipation  coefficient  needed  to  be  set  to  a  value  of  1.0  to  stabilize  the  calculation  in 
the  presence  of  the  very  strong  primary  shock  which,  at  early  times,  had  a  Mach  number  of 
~  15. 

The  computational  grids  were  assembled  from  two  zones  in  which  the  axial  and  radial 
node  distributions  were  generated  independently  using  one  of  Roberts’  [28]  stretching  trans¬ 
formations.  (See  also  [29],  Section  5-b.l).  In  each  zone  there  were  N  -f  1  nodes,  including 
the  end  points,  located  at 

z  =  zaf}  +  zb(l  -  rj)  ,  (8) 

where 

[(/3  +  2a)A  —  J3  +  2a] 

V  =  (2a  -f  1)  [1  +  A]  ’ 
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Axial  Zone 

a 

P 

N 

Wall  Condition 

a 

P 

N 

0  <  x  <  0.2m 

0 

1.04 

100 

(b)  Inviscid 

0.5 

100 

50 

0.2  <  x  <  2.187m 

0 

10.0 

494 

Viscous 

0.5 

1.1 

80 

Table  4:  Grid  generation  parameters,  (a)  Axial  distribution,  (b)  Radial  distribution. 

J^l) 

V  =  J  -  0 ...N  . 

Details  of  the  zonal  boundaries  and  stretching  parameters  are  given  in  Table  4.  The  same 
axial  distribution  was  used  for  both  inviscid-wall  and  viscous-wall  solutions.  The  axial 
resolution  of  the  grids  was  determined  by  the  need  to  ensure  conservation  in  the  initial 
conical  expansion  and  then  to  retain  adequate  resolution  of  the  starting  shock  structure  as 
the  structure  was  swept  through  the  nozzle.  The  radial  grid  distribution  for  the  viscous-wall 
solution  was  quite  different  from  that  of  the  inviscid  solution.  While  a  uniform  distribution 
of  N  =  50  points  was  found  to  be  adequate  for  the  inviscid  solution,  the  viscous  grid 
had  N  =  80  and  was  stretched  to  cluster  nodes  close  to  the  wall  to  resolve  the  boundary 
layer.  However,  the  number  of  nodes  used  and  the  clustering  specified  was  something  of  a 
compromise  and  the  calculation  did  not  adequately  resolve  the  boundary  layer  in  the  early 
part  of  the  expansion. 

It  should  also  be  noted  that  (1)  zero-gradient  boundary  conditions  were  implemented 
as  zero-order  extrapolations  and  (2)  zero-gradient  boundary  conditions  at  the  wall  were 
implemented  by  specifying  the  derivatives  with  respect  to  the  y-coordinate  (rather  than 
normal  to  the  wall).  Future  calculations  will  be  performed  on  a  modified  grid  with  grid  lines 
normal  to  the  wall. 

7.1  Results 

Figures  19  and  20  show  sequences  of  Mach  number  contours  for  the  inviscid-wall  and  viscous- 
wall  simulations,  respectively.  The  initial  pressure  Pin;t  =  150Po  is  higher  than  typically 
used  in  the  facility,  but  was  the  smallest  value  that  could  be  used  with  moderate  levels  of 
artificial  diffusion.  As  a  result,  the  unsteady  expansion  was  only  evident  in  the  solutions  up 
to  t  =  0.3ms.  See,  for  example,  the  ’’turning  up”  of  the  Mach  contours  at  t  =  0.3ms  and 
x  —  0.9m  in  Fig.  19.  Beyond  that  time,  the  upstream-facing  shock  became  the  limiting 
wave  for  steady  flow. 
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The  most  obvious  multidimensional  effect  shown  in  the  inviscid-wall  solution  is  the  dis¬ 
tortion  of  the  contact  discontinuity  and  the  upstream-facing  shock  into  waves  that  are  not 
planar.  The  distortion  is  enough  to  delay  the  flow  settling  time  by  approximately  0.04ms  be¬ 
yond  that  computed  in  the  Q1D  results.  This  is  approximately  the  size  of  the  delay  between 
the  Q1D  results  and  the  experimental  values  as  shown  in  Fig.  17.  Although  the  contact 
discontinuity  was  distorted  severely,  the  primary  shock  remained  nearly  planar  and  its  speed 
matched  that  computed  in  the  Q1D  results  (cf.  Fig.  15). 

The  introduction  of  viscous  effects  caused  a  thick  boundary  layer  to  develop  along  the 
nozzle  wall  (see  Fig.  20).  The  starting  shock  structure  remained  essentially  the  same, 
with  the  primary  shock  again  nearly  planar.  At  t  =  0.511ms  the  primary  shock  arrived 
at  x  =  1.8m  and  had  a  speed  of  3214 m/s.  The  upstream-facing  shock  arrived  later  at 
t  =  0.643ms  with  a  speed  of  2643 m/s.  The  low  pressure  gas  which  passed  through  the 
primary  shock  and  accumulated  before  the  contact  discontinuity,  had  a  Mach  number  of 
approximately  2.  This  is  fairly  close  to  the  one- dimensional  Mach  number  limit  for  an 
infinitely  strong  shock  compression  (see  e.g.  [30]) 

(  2  V/2 

M=fcj)  -1'89’ 

with  7  =  1.4. 

Further  details  of  the  starting  shock  structure  are  shown  in  figures  21  and  22.  The  static 
pressure  contours  identify  the  primary  shock  and  the  upstream-facing  shock.  Alternatively, 
the  temperature  field  identifies  the  contact  discontinuity  and  the  boundary  layer  with  its 
large  temperature  gradients.  The  static  pressure  is  essentially  constant  across  both  of  these 
features.  The  density  field,  being  sensitive  to  both  pressure  and  temperature,  is  more  com¬ 
plicated  and  shows  all  of  the  flow  features.  The  velocity  field  shows  abrupt  changes  at  the 
two  shocks  and  complex  recirculation  regions  between  them.  It  also  shows  the  formation  of 
the  boundary  layer  (for  a:  <  1.75  at  t  =  0.6ms)  by  the  entrainment  of  gas  along  the  wall  [31] 
(23). 

Figure  23  shows  the  evolution  of  the  total  boundary-layer  thickness  6  (as  measured  by 
the  99%  point  in  the  total  enthalpy  H  —  h  -f  0.5u2  profile),  the  displacement  thickness  6* 
and  the  momentum  thickness  6  as  defined  by  the  expressions  [32] 


r%  ,  frW<xll  ( 

rwall  ~  (ru>all  ~  ^  )  =  2  /  I  1 


r  wall 


-(%--*}■  =  2  -CL 
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PU)rdr> 

(10) 
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The  total  thickness  and  the  displacement  thickness  settle  to  a  steady  state  shortly  after  the 
passage  of  the  upstream-facing  shock.  The  momentum  thickness,  however,  takes  a  signifi¬ 
cantly  longer  time  to  settle.  At  t  —  0.6ms,  the  upstream-facing  shock  has  reached  x  ~  1.7m 
while  only  the  first  0.6m  of  the  momentum  thickness  profile  has  reached  steady  state.  This 
agrees  reasonably  well  with  the  settling  time  correlations  reported  in  [33]. 

By  t  =  1.2ms  the  boundary  layer  has  approached  steady  state  and  8  appears  to  be  a 
linear  function  of  x,  reaching  8  ~  46mm  at  x  =  1.8m.  The  displacement  thickness  remains 
small  for  x  <  1.0m  but  then  grows  linearly  (in  x )  to  reach  approximately  0.5£  at  the  nozzle 
exit.  The  small  initial  growth  rate  for  x  <  0.8m  is  due  to  the  combination  of  the  hot  free- 
stream  conditions  and  the  relatively  cold  wall  (300° A').  Momentum  thickness  grows  most 
rapidly  while  the  displacement  thickness  is  small  but  levels  off  once  8*  begins  to  grow  rapidly. 
At  x  —  1.8m  the  displacement  thickness  and  momentum  thickness  are  approximately  23mm 
and  3.2mm,  respectively. 

The  Pitot  pressure,  as  computed  by  Eqn.(7),  is  shown  in  Fig.  24.  Although  the  traces  are 
somewhat  noisy,  the  starting  shock  structure  is  obvious.  The  settling  of  the  traces  (r  =  112) 
from  the  passage  of  the  upstream- facing  shock  to  t  =  0.8ms  is  probably  caused  by  the  set¬ 
tling  of  the  boundary  layer  along  the  nozzle  wall.  Unfortunately,  the  large  and  persistent 
disturbances  seen  in  the  experimental  data  are  not  seen  even  in  these  two-dimensional  sim¬ 
ulations.  Unsteady  inflow  conditions  and  fully  three-dimensional  flow  models  have  yet  to  be 
investigated. 

The  radial  distributions  of  Pitot  pressure  are  compared  with  quasi-steady  experimen¬ 
tal  values  in  Fig.  25.  At  x  =  1.8m  the  computed  values  for  the  viscous  wall  simulation 
fortuitously  capture  the  spurious  data  points  at  the  nozzle  axis.  Although  the  bulk  of  the 
computed  values  are  reliable  and  match  the  experimental  values  well,  those  near  the  nozzle 
axis  are  subject  to  a  particularly  noisy  boundary  condition  and  should  not  be  trusted.  At 
x  =  2.0m  there  is  obvious  disagreement  with  the  experiment  near  the  axis.  The  experimen¬ 
tal  results  also  indicate  that  the  physical  boundary  layer  is  thicker  than  that  computed  in 
the  viscous-wall  simulation.  Note  that  the  simulation  was  for  laminar  flow  only  whereas  the 
physical  flow  may  have  been  partly  turbulent. 

8  CONCLUDING  REMARKS 

From  the  experimental  point  of  view,  the  nozzle  appears  to  produce  a  satisfactory  test  flow 
once  it  has  settled  to  a  (quasi-)steady  state.  Hence,  the  simple  inviscid  design  approach  is 
adequate  but  would  benefit  from  the  addition  of  a  boundary  layer  correction.  It  may  be 
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possible  to  partially  compensate  for  this  oversight  by  manufacturing  throats  with  smaller 
diameters  while  maintaining  the  current  dimensions  of  the  contoured  expansion. 

Although  the  steady  performance  is  adequate,  the  settling  times  are  much  longer  than 
predicted  by  the  classical  quasi-one-dimensional  models.  This  delay  in  flow  establishment  is 
especially  important  at  high  enthalpy  ( H ,  >  16 MJ/kg)  operating  conditions  as  the  nozzle 
supply  conditions  arc  maintained  for  relatively  short  times  after  shock  reflection.  Fortunately, 
this  problem  is  compensated  by  the  magnitude  of  the  disturbances  becoming  small  at  higher 

H.. 

From  the  computational  point  of  view,  the  nozzle  starting  processes  are  not  fully  un¬ 
derstood.  The  quasi-one-dimensional  model  appears  to  provide  an  adequate  description  of 
the  starting  shock  structure  but  fails  to  predict  the  long  settling  times  observed  in  the  ex¬ 
periment.  The  axisymmetric  simulations  show  that,  while  the  primary  shock  remain  nearly 
planar,  the  contact  discontinuity  and  the  upstream-facing  shock  are  distorted.  This  distor¬ 
tion  is  responsible  for  slightly  delaying  their  arrival  at  the  nozzle  exit  plane  (beyond  the  times 
indicated  by  the  quasi-one-dimensional  model).  However,  the  axisymmetric  simulations  did 
not  provide  any  further  insight  to  the  mechanism  behind  the  persistent  disturbances  seen 
in  the  experimental  data.  Further  simulations  need  to  be  performed  which  may  include  (a) 
the  dynamics  of  the  nozzle  supply  region  (i.e.  the  shock  reflection  region  of  the  shock  tube); 
and/or  (b)  fully  three-dimensional  effects  (starting  with  a  slightly  elliptic  throat). 
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Notation 

a  :  speed  of  sound 

A  :  interface  area  (Q1D  code) 

c  :  thermodynamic  polynomial  coefficients  (Appendix  A) 

d  :  diameter 

e  :  specific  internal  energy 

E  :  total  energy 

/  :  density  weighting  parameter  (QlD  code) 

H  :  total  enthalpy 

L  :  nozzle  length 

M  :  Mach  number 

N  :  number  of  computational  nodes/cells 

P  :  pressure  . 

r  :  radial  coordinate 

T  :  temperature 

t  :  time 

u  :  axial  velocity 

U  :  Riemann  invariant  (QlD  code) 

v  :  radial  velocity 

V  :  cell  volume  (QlD  code) 

w  :  wave  velocity  (QlD  code) 

x  :  axial  coordinate 

X,  Y  :  thermodynamic  properties  (Appendix  A) 

Z  :  intermediate  quantity  (QlD  code) 

a  :  transient  inflow  parameter 

a,  :  coordinate  stretching  parameters 

7  :  ratio  of  specific  heats 

St  :  time  delay  for  normalizing 

8  :  boundary  layer  thickness 

8*  :  boundary  layer  displacement  thickness 

A  :  intermediate  coordinate;  driver  tube  compression  ratio 
p  :  density 

t  :  normalized  time  (transient  inflow  parameter) 

rj  :  radial  computational  coordinate 

£  :  axial  computational  coordinate 
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Subscripts: 

* 

CL,  0 j  C|  d 

A,  B,  C ,  D 

d 

e 

L,R 

3 

s  -  ref 
t 

wall  . •• 


throat  condition 

labels  (transient  inflow  description) 

locations  in  the  inviscid  design 

inviscid  design  value 

free-stream  conditions 

left  and  right  states  (Q1D  code) 

nozzle  supply/stagnation  condition 

reference  supply/stagnation  condition 

truncated  nozzle  value 

value  at  the  nozzle  wall 


Superscripts: 

*  :  intermediate  states  (Q1D  code) 
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Y 

Units 

X 

Units 

Y 

vvmm 

v 

max 

Extrapolation 

hje 

- 

e 

MJ/kg 

0.5223 

4.7609 

constant 

T 

K/ 1000 

e 

MJ/kg 

0.5223 

4.7609 

linear 

e 

MJ/kg 

T 

K/im 

0.7215 

5.0 

linear 

Y 

X 

c0 

Cl 

c2 

h/e 

e 

0.441746 

-1.014804  x  lO"1 

3.458821  x  lO"2 

T 

e 

3.370113  x  10-2 

1.426067 

-0.236355 

e 

T 

4.937905  x  10"3 

0.614562 

0.164100 

Y 

X 

c3 

c4 

c5 

h/e 

e 

-6.077331  x  10"3 

4.211913  x  10"4 

0 

T 

e 

6.481913  x  lO"2 

-9.236112  x  10"3 

5.254391  x  10~4 

e 

T 

-3.273006  x  10-2 

3.273170  x  lO"3 

-1.196175  x  lO"4 

Table  5:  Fitted  ranges  and  coefficients  for  the  thermodynamic  properties. 

A  Thermodynamic  Model 

The  thermodynamic  model  used  in  the  quasi-one-dimensional  code  assumed  a  non-reacting 
mixture  of  nitrogen  and  oxygen  with  mass  fractions  0.7686  and  0.2314,  respectively.  The 
specific  heats  for  each  component  were  assumed  to  be  functions  of  temperature  only,  and 
the  curve  fits  in  [22]  were  used  to  provide  a  small  data  base  for  the  gas  mixture.  Polynomial 
expressions  for  h/e  —  P/(pe)  +  1,  temperature,  and  internal  energy  were  fitted  to  this  data 
base.  The  generic  polynomial  was 

V  =  co  +  C\X  +  C2-X2  +  C3A  3  -f  C4X4  +  C5 X5  ,  (12) 

and  the  fitted  coefficients  are  given  in  Table  5.  Outside  their  fitted  ranges,  estimates  oi  h/e 
were  held  constant  while  T  and  e  were  linearly  extrapolated. 
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B  Quasi- One-Dimensional  Formulation 


The  governing  equations  for  quasi-one-dimensional  inviscid  flow  without  heat  addition  are 
formulated  with  respect  to  the  finite  volume  shown  in  Fig.  26.  The  geometry  is  specified  as 
a  set  of  N  cells  separated  by  interfaces  located  at  Xj+ 1 ,  j  =  0...N  with  areas  Aj+i.  Note 
that  Xj+i  =  x(j+1)_i .  The  volume  of  cell  j  is 

Vi  =  5  +  (13) 


The  flow  solution  is  recorded  as  cell  averages  of  density  pj,  velocity  Uj,  and  specific  internal 
energy  ey 

The  differential  equations  governing  the  conservation  of  mass,  momentum,  and  total 
energy  are  approximated  by 


(14) 


+pi(Aj+ 1  -  Aj_i) 


(15) 


\(P,+k  +Pi+iEi+  0ui+lAj+k 


(16) 


where  the  total  energy  E  =  e  +  ^u2,  and  P  =  P(p,  e)  is  the  equation  of  state  described  in 
Appendix  A.  These  conservation  equations  are  integrated  forward  in  time  with  a  predictor- 
corrector  technique,  while  interface  values  (with  subscripts  ±^)  are  evaluated  in  a  stable 
manner  with  an  Osher-type  Riemann  solver  [25]. 

At  the  start  of  each  time  step,  the  state  of  the  flow  (consisting  of  a  set  of  values  for 
p,  u ,  e,  P,  a )  either  side  of  each  interface  is  interpolated  from  the  set  of  cell  averaged  states 
by  assuming  piecewise  linear  variation  of  the  variables  within  cells.  The  slopes  of  these 
linear  segments  are  limited  by  the  nonlinear  SUPER  REE  limiter  [26].  The  Riemann  solver 
then  assumes  that  a  spatially  constant  left  state  (subscript  L)  and  right  state  (subscript 
R)  interact  through  a  pair  of  finite-amplitude  compression  or  rarefaction  waves.  Perfect  gas 
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relations  [34]  are  used  to  obtain  the  states  ( L *,  R*)  in  the  gas  after  the  passage  of  left-moving 
and  right-moving  waves,  respectively.  The  expressions  implemented  in  the  code  are 

PL-P«~  Pl  2at(l  +  Z)  ]  ’  (17) 


Ur  —  U  d  — 


ULZ  +  Ur 
1  +  Z  ’ 


where  the  Riemann  invariants  are 


Ul  =  ul+  ,  and 

7  -  1 

7 7  2afl 

Ur  =  uR - -  , 

7  -  1 

and  the  intermediate  variable  Z  is  given  by 


O-L  \Pr) 


Given  the  pressure  and  velocity  in  the  regions  behind  the  waves,  the  other  flow  properties 
may  be  evaluated.  Sound  speeds  are  evaluated  from  the  invariants  (19)  as 

*1  =  (Ul  ~  ul)  (7  -  l)/2  ,  and 

a%  =  K-F«)(7-l)/2-  (21) 

The  specific  internal  energy  is  obtained  from 

.  - 

L  ~  +T^T)' 


el  = 


M2 

7(7  - 1)  ’ 


and  the  density  is  obtained  from  the  equation  of  state  as 


.  Pr 


If  the  pressure  rises  across  the  left-moving  wave,  it  is  assumed  to  be  a  shock  and  its 
velocity  relative  to  the  initial  left  state  is  given  by 


uL-wL  = 


7+1  Pl  (PI  7-1 
2  pL\PL  7  +  I 
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else,  it  is  assumed  to  be  a  rarefaction  with  relative  velocity 


ul  -  wl  =  aL  .  (25) 

Here  wl  is  the  velocity  of  the  wave  relative  to  the  cell  boundaries.  A  similar  procedure  is 
used  to  obtain  the  velocity  of  the  right-moving  wave. 

In  the  preceding  perfect  gas  equations  (17)  -  (24),  an  effective  7  was  used  to  include  real 
gas  effects  in  an  approximate  manner  and  was  evaluated  as  a  density  weighted  function  [35] 


where 


and 


7  =  flL  +  (1  -  /)7K . 


/  = 


\fpi 

\/PL  +  y/]>R  ’ 


7  L 
lR 


Pl 

pL  &L 


+  1 


Pr 

PRCR 


+  1  ■ 


(26) 

(27) 


(28) 


Given  the  signs  and  relative  magnitudes  of  the  wave  speeds,  the  interface  conditions  used 
in  the  conservation  equations  (14)  -  (16)  may  be  selected  or  interpolated  from  the  4  flow 
states  using  the  following  logic. 
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if  ( u *  >  0)  then 


The  contact-discontinuity  has  moved  to  the  right 
and  the  interface  state  is  determined  from  the 
L  and  L*  states. 

if  (P*  >  Pi)  then 
The  left-moving  wave  is  a  shock, 
if  (wi  >  0)  then 

All  waves  have  moved  to  the  right. 

Interface  state  is  equal  to  L. 
else 

Interface  state  is  equal  to  L* . 
endif 
else 

The  left-moving  wave  is  a  rarefaction, 
if  ( ui  —  ai  >  0)  then 
All  waves  have  moved  to  the  right. 

Interface  state  is  equal  to  L. 
elseif  ( u*L  —  a*L  >  0)  then 
The  rarefaction  straddles  the  interface. 
Interpolate  the  interface  state  from 
states  L  and  L* . 
else 

The  entire  rarefaction  moved  to  the 
left  of  the  interface. 

Interface  state  is  equal  to  L*. 
endif 
endif 


The  contact  discontinuity  has  moved  to  the  left 
and  the  interface  state  is  determined  from  the 
R  and  R*  states  in  a  similar  manner... 

endif 
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Figure  1:  Reflected-shock  tunnel  operation,  (a)  Tunnel  schematic,  (b)  x  —  t  wave  diagram. 
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Figure  2:  Typical  history  of  the  (unfiltered)  nozzle  supply  pressure  showing  the  principal 
events:  (a)  arrival  of  the  incident  shock;  (b)  reflected  shock;  (c)  establishment  of  equilibrium 
pressure;  (d)  driver  gas  contamination. 
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Figure  3:  Quasi-one-dimensional  model  of  the  nozzle  starting  process  in  both  the  physical 
plane  and  the  x  —  t  plane.  Labelled  features  are:  (a)  primary  shock;  (b)  contact  surface;  (c) 
upstream-facing  shock;  (d)  upstream  head  of  the  unsteady  expansion;  (e)  steady  expansion. 
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Figure  4:  Schematic  view  of  a  contoured  nozzle  showing  regions  of:  (1)  stagnant  nozzle 
supply  gas;  (2)  transition  to  source  flow;  (3)  source  flow;  (4)  transition  to  parallel  flow;  (5) 
parallel  flow. 


Figure  5:  Characteristic  mesh  for  regions  3,  4  and  5.  The  interpolated  streamline  extends 
from  A  to  D.  Note  that  the  axial  coordinate  has  been  compressed. 
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Figure  8:  Filtered  Pitot  pressure  histories  at  x  ~  1.8m  for  shots  1405  and  1406.  Note  that 
the  time  base  has  been  shifted  such  that  shock  reflection  occurs  at  t  =  0. 


V/277J,  km/s 


Figure  10:  Times  at  which  events  occur  at  the  nozzle  exit  plane  for  a  4mm  primary  diaphragm 
and  A  =  60.  A.  primary  shock  arrival  (Q1D  model);  B.  settling  time  (10%,  Q1D  model); 
C.  contamination  by  driver  gas;  x  primary  shock  arrival  (experimental);  O  settling  time 
(experimental);  O  passage  of  the  expansion  fan  (experimental). 
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Figure  12:  Flow  geometry  for  the  numerical  simulations. 
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Figure  13:  Transient  supply  pressure, 
transient  inflow  approximation. 


Dashed  line:  filtered  Pt  for  shot  1405;  Solid  line: 
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(a)  Pinit  =  33 Pa  (b)  Pinii  =  150 Pa 


0.0  0.5  1.0  1.5  2.0  0.0  0.5  1.0  1.5  2.0 

x,  m  x,  m 


Figure  15:  Evolution  of  the  density  profile  for  (a)  Pinit  =  33Pa,  (b)  Pinit  =  150Pa.  Solid 
lines  are  for  the  times  as  shown.  Dashed  lines  show  the  solution  at  t  ~  1.2ms. 
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Figure  16:  Wave  diagrams  for  various  initial  pressures  showing:  (a)  primary  shock;  (b) 
contact  surface;  (c)  upstream-facing  shock;  (d)  upstream  head  of  the  unsteady  expansion. 


Figure  17:  The  variation  of  event  times  with  initial  pressure.  O  settling  time  (10%)  for 
steady  inflow;  *  settling  time  (10%)  -  transient  inflow;  Q  primary  shock  passage  for  steady 
inflow;  A  primary  shock  passage  for  transient  inflow.  Hatched  areas  show  the  estimated 
experimental  times  for  (a)  primary  shock  arrival  and  (b)  arrival  of  the  upstream-head  of  the 
unsteady  expansion. 
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Figure  18:  Comparison  of  (a)  raw  and  (b)  normalized  Pitot  pressure  histories  and  experi¬ 
mental  data  for  r  =  7.5 mm,  shot  1406:  solid  line  =  steady  throat  conditions;  short  dash  = 
transient  throat  conditions;  long  dash  =  experiment. 
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Figure  19:  Mach  number  contours  showing  the  axisymmetric  starting  process  with  an 
inviscid-wall  boundary  condition.  (Contour  levels  every  0.5.) 
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Figure  20:  Mach  number  contours  showing  the  axisymmetric  starting  process  with  a  viscous- 
wall  boundary  condition.  (Contour  levels  every  0.5.) 
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Figure  22:  Details  of  the  shock  structure  for  the  viscous-wall  simulation  at  t  =  0.6ms.  (a) 
Mach  number,  contour  levels  every  0.5;  (b)  velocity  field  in  the  frame  of  reference  stationary 
with  respect  to  the  primary  shock  ( u  ~  3214 m/s). 
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Figure  23:  Boundary  layer  development  for  the  viscous  wall  simulation:  solid  line  =  6\  long 
dash  =  8 *;  short  dash  =  9. 
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Figure  24:  Comparison  of  normalized  Pitot  pressure  histories  for  the  viscous-wall  simulation 
and  shot  1406:  solid  line  =  simulation;  dash  =  experiment. 


Figure  25:  Quasi-steady  normalized  Pitot  pressures  across  the  nozzle  exit  plane  ( x  ~  1.8m) 
and  further  downstream  at  x  ~  2.0m:  solid  line  =  viscous-wall  simulation;  dash  =  inviscid- 
wall  simulation. 


Figure  26:  The  generic  cell  used  in  the  quasi-one-dimensional  formulation. 
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